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I.  INTRODUCTION 


Accurate  and  unbiased  surface  wave  measurements  are  required  if 
surface  waves  are  to  be  used  for  treaty  monitoring  purposes. 
Conventional  time  domain  Ms  measurements  provide  a  rough  estimate  of 
the  size  of  the  seismic  source.  Differences  in  travel  path,  however, 
cause  variations  in  Mg  that  are  unacceptable  in  cases  where  0.1 
magnitude  unit  may  be  important.  On  short  paths  (less  than  about  30 
degrees),  the  seismogram  is  dominated  by  the  Airy  phase  and  20  second 
amplitudes  are  difficult  or  impossible  to  measure  accurately.  On  long 
paths,  the  surface  waves  are  we  1 1 -dispersed  but  attenuation 
differences  can  cause  large  amplitude  variations.  Since  explosion 
yield  estimation  is  usually  done  by  defining  a  magnitude  versus  yield 
curve  at  one  or  more  test  sites,  path  differences  and  differences  in 
station  coverage  can  lead  to  a  biased  result  when  using  this  curve  to 
estimate  yields  at  an  uncalibrated  test  site. 

In  this  report,  we  describe  a  method  for  obtaining  more  stable 
and  reliable  surface  wave  measurements  by  using  all  of  the  information 
contained  in  the  surface  wave  to  estimate  the  Green’s  function  for  the 
path.  We  use  Rayleigh  waves  from  explosions  at  the  Nevada  Test  Site 
(NTS)  and  the  Soviet  East  Kazakh  test  site  to  find  Green’s  functions 
for  the  paths  between  these  test  sites  and  SRO  and  WWSSN  stations. 
The  phase  of  the  Green’s  function  can  be  used  to  construct  a  phase- 
matched  filter,  which  may  be  used  to  recover  spectral  amplitudes  for 
surface  waves  from  other  events  observed  along  nearby  paths.  The 
moment  of  the  explosion  can  then  be  found  by  dividing  the  spectral 
amplitude  by  the  amplitude  of  the  Green’s  function  for  the  path. 

Figure  1  illustrates  the  difficulty  in  measuring  surface  wave 
amplitudes  using  standard  methods.  The  top  figure  is  a  Rayleigh  wave 
from  an  NTS  explosion  recorded  at  Station  GOL,  at  a  distance  of  977 
km.  At  this  close  distance  many  frequencies  arrive  at  the  same  time, 
so  that  the  seismogram  is  dominated  by  the  Airy  phase.  A  reliable  20 
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Seismogram  and  spectrum,  NTS-GOL 


This  figure  illustrates  the  difficulty  in  measuring  surface  wave 
amplitudes  using  standard  methods.  The  top  figure  is  a  Rayleigh  wave 
from  the  NTS  explosion  SANDREEF  recorded  at  Station  GOL.  The 
seismogram  is  dominated  by  the  Airy  phase,  so  a  reliable  20  second 
amplitude  cannot  be  measured.  The  lower  figure  shows  the  average 
spectrum  (solid  line)  and  *  1  standard  deviation  (dotted  line)  for  24 
spectra  at  Station  GOL.  The  amplitude  depends  strongly  on  the 
frequency  measured. 


second  amplitude  cannot  be  measured.  Measuring  the  spectral  amplitude 
instead  of  the  seismogram  removes  this  problem.  The  lower  figure 
shows  the  average  spectrum  (solid  line)  and  +  1  standard  deviation 
(dotted  line)  for  24  spectra  at  Station  COL.  Another  problem  still 
remains,  however.  The  amplitude  depends  strongly  on  the  frequency 
measured,  changing  by  a  factor  of  2  between  0.035  Hz  and  0.05  Hz,  and 
by  a  factor  of  5  between  0.02  Hz  and  0.05  Hz.  Figure  2  shows  the 
amplitude  and  phase  of  the  spectrum  after  the  path  correction  has  been 
applied.  Notice  that  in  addition  to  correcting  for  path  effects,  the 
path  correction  is  also  a  frequency  correction,  flattening  the 
spectrum  over  the  frequency  band  of  interest. 

The  procedures  recommended  in  this  report  are  based  on  our 
experience  applying  surface  wave  analysis  methods  to  real  and 
synthetic  data  over  the  past  few  years.  Bache,  et  al.  (1978a)  used 
inversion  of  surface  wave  dispersion  to  estimate  the  path  structure 
between  NTS  and  WWSSN  Stations  ALQ  and  TUC.  Bache,  et  al.  (1978b) 
used  these  structures  to  estimate  the  source  size  of  NTS  explosions. 
Rodi,  et  al.  (1978)  inverted  Rayleigh  waves  from  explosions  at  the 
French  Sahara  test  site  to  estimate  path  structure  in  Africa.  Our 
procedure  is  similar  to  that  used  by  Bache,  et  al.  (1978a),  except 
that  we  use  phase-matched  filtering  to  recover  phase  velocities  and 
spectral  amplitudes;  we  invert  for  the  Q  structure  of  the  path  as  well 
as  the  elastic  structure;  and  the  programs  have  been  streamlined  so 
many  paths  can  be  processed  in  an  efficient  manner.  Tests  of  the 
phase-matched  filtering  and  inversion  programs  are  given  by  Wang,  et 
al.  (1981).  The  computer  programs  used  in  this  study  are  described  in 
detail  in  a  report  by  Stevens,  et  al.  (1982). 
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Corrected  Amplitude  6f~  Corrected  Phase 


Figure  2.  Path  corrected  amplitude  and  phase  for  the  path  NTS-GOL.  The  variation  in  amplitude  with 
frequency  is  greatly  reduced.  The  phase  correction  can  be  used  to  determine  phase  shifts 
(such  as  Rayleigh  wave  reversals)  between  events. 


II.  RAYLEIGH  WAVES  FROM  AN  EXPLOSION  SOURCE 


Our  objective  is  to  find  the  Green’s  function  for  a  given  source 
to  receiver  path  which  can  be  used  to  determine  the  explosion  moment. 

Our  approach  to  the  problem  has  been  to  assume  that  the  average  path  ' 
structure  can  be  adequately  modeled  by  a  plane- layered,  attenuating 
elastic  medium.  Since  the  explosions  are  located  in  a  common  source 
region,  but  the  Rayleigh  waves  travel  along  different  paths,  we  also 
use  an  approximation  which  allows  us  to  separate  the  spectrum  into 
source  region  and  path  dependent  functions  (Bache,  et  aL.,  1978a). 
This  approximation  assumes  energy  conservation  and  no  mode  conversion 
when  a  Rayleigh  wave  crosses  a  boundary  between  adjacent  p I ane- I ayered 
structures.  The  accuracy  of  this  approximation  has  been  demonstrated 
for  selected  cases  by  comparison  with  finite  difference  calculations 
by  Glover  and  Harkrider  (1981).  We  also  assume  an  isotropic  explosion 
source,  neglecting  non  isotropic  source  effects  such  as  explosion- 
induced  block  motion  and  tectonic  strain  release.  In  the  construction 
of  the  path  Green’s  functions,  we  try  to  use  events  with  little 
evidence  of  non  isotropic  source  effects. 

Given  these  assumptions,  we  can  write  the  vertical  component  of 
the  Rayleigh  wave  as  a  product  of  functions  each  of  which  depends  only 
on  the  source,  the  source  region,  or  the  path: 

u(w,r)  ,  M0'S1(w)S2(w)exp(-72r)/(^7R^  exp(<  (^-wr/^)]  (1) 

where  Mq'  is  the  normalized  moment  of  the  explosion  discussed  below, 
r  is  the  source  to  receiver  distance,  ae  is  the  radius  of  the  earth, 

C2  is  the  frequency  dependent  phase  velocity  for  the  path  structure, 

72  is  the  frequency  dependent  attenuation  coefficient  due  to  the 
anelastic  structure  of  the  path,  and  Iq  is  the  initial  phase  which  is 
equal  to  -3*/4  for  an  explosion  with  step  function  moment.  Sj(w)  and 
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$2(u)  are  real  functions  that  depend  on  the  earth  structure  of  the 
source  region  and  path  respectively.  S^(v)  is  positive  for  a  shallow 
source.  Using  the  notation  of  Harkrider  (1964,  1970), 


$2  (u)  ~  c2  ^  (2) 

E|  and  E3  are  the  eigenfunctions  of  horizontal  displacement  and  normal 
traction  respectively,  and  are  written  in  Harkrider's  notation  as  E^  = 
Us(h)/£o,  E3  =  a*(h)/(£g/c)  is  the  shear  modulus  at  the  source). 
They  are  functions  of  the  source  depth  h  and  the  large  scale 
properties  of  the  earth  structure,  c^  and  C2  are  the  phase  velocities 
of  the  source  and  path  structures,  and  Arj  and  Ar2  are  the  frequency 
dependent  Rayleigh  wave  amplitude  functions  for  the  source  and  path 
regions  which  depend  on  the  large  scale  properties  of  the  earth 
structure. 

In  the  following  analysis,  the  source  region  elastic  structure 
is  assumed  to  be  constrained  by  independent  geophysical  evidence  and 
the  same  structure  will  be  used  for  all  events  at  a  given  test  site. 
The  quantities  which  must  be  determined  in  Equation  (1)  are  C2,  S2, 
and  72. 

The  normalized  moment  Mg  is  equal  to  3/J^/a^Mg,  where  Mg  is  the 
conventional  definition  of  the  moment  and  a  and  p  are  the 
compress i ona I  velocity  and  shear  velocity  at  the  source.  For  a 
Poisson  solid  Mg  a  Mg.  We  define  Mg  this  way  so  that  it  will  be  a 
measure  of  the  strength  of  surface  wave  excitation  by  the  source,  and 
so  the  Green's  function  will  be  independent  of  local  material 
properties.  With  this  definition,  two  nearby  explosions  in  media  with 
different  Poisson's  ratios  will  produce  surface  waves  of  the  same 
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III.  EXPLOSION  DATA  FOR  THE  NEVADA  TEST  SITE  AND  EAST  KAZAKH 


The  data  used  in  this  study  were  Rayleigh  wave  seismograms  from 
underground  nuclear  explosions  at  the  Nevada  Test  Site  (NTS)  and  the 
Soviet  East  Kazakh  test  site  (EKZ) .  Stations  were  selected  which  had 
good  signa l-to-noise  ratios  for  at  least  two  events.  The  location  of 
the  recording  stations  relative  to  the  test  sites  are  shown  in  Figures 
3,  4,  and  5.  Twenty-two  explosions  at  the  East  Kazakh  test  site  were 
recorded  at  12  SRO  stations.  Forty-eight  NTS  explosions  were  recorded 
at  12  SRO  stations  and  24  WWSSN  stations. 

The  WWSSN  data  were  hand  digitized  by  ENSCO  Incorporated,  and 
installed  at  the  Center  for  Seismic  Studies  in  Arlington,  Virginia. 
Except  for  occasional  timing  errors,  we  found  the  WWSSN  data  to  be  of 
very  good  quality,  and  we  have  had  little  difficultly  in  applying  the 
inversion  procedure  to  this  data.  All  of  the  data  processing  was  done 
on  the  VAX  11/780  at  the  Center  for  Seismic  Studies. 

The  paths  shown  in  Figures  3  through  5  cover  a  variety  of  earth 
structures  and  a  considerable  range  of  distances.  The  paths  from  NTS 
to  WWSSN  stations  are  relatively  short,  ranging  from  900  to  4400  km 
and  are  all  continental  paths.  The  SRO  stations  are  all  quite  distant 
from  NTS,  ranging  from  5600  to  14,000  km.  Three  of  these  paths  are 
almost  entirely  oceanic,  while  the  remainder  are  mixed  continental  and 
oceanic.  Most  of  the  paths  from  East  Kazakh  to  SRO  stations  are 
between  2000  and  5000  km,  but  two  of  the  paths  are  longer  than  10,000 
km.  Seismograms  are  nearly  identical  for  all  events  observed  along 
the  same  path,  except  for  anomalous  events  with  reversed  Rayleigh 
waves . 
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Figure  4.  Surface  wave  paths  from  NTS  to  12  SRQ  stations.  Stations  farther 
than  90  degrees  from  the  test  site  are  shown  on  the  outside  of  the 
figure  in  the  direction  of  the  station. 
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IV.  DATA  ANALYSIS  AND  INVERSION  PROCEDURE 


A  block  diagram  of  the  method  used  for  the  data  processing  is 
shown  in  Figure  6.  We  start  by  finding  a  set  of  Rayleigh  wave 
seismograms  from  a  given  source  region  and  recorded  at  the  same 
station  which  show  consistent  signals  and  have  good  s igna I -to-noi se 
ratios.  We  then  find  the  phase  velocities,  group  velocities  and 
amplitude  spectra  for  each  of  these  seismograms  using  narrow  band 
filtering  (Archambeau,  et  at.,  1966),  Dziewonski,  et  at.,  1969, 
Savino,  et  at.,  1980)  and  phase-matched  filtering  (Herrin  and  Goforth, 
1977;  Wang,  et  at.,  1981).  Group  velocities  are  found  first  by 
applying  a  set  of  narrow  band  filters  at  successive  center  frequencies 
and  computing  the  arrival  times  of  the  energy  envelopes.  The  phase 
velocities  are  found  next  by  integrating  the  group  delays  to  find  the 
phase  as  a  function  of  frequency,  forming  a  phase-matched  filter, 
cross-correlating  the  filter  with  the  seismograms,  windowing  the 
compressed  signal,  and  unwrapping  the  residual  phase  spectrum  to 
correct  the  initial  phase  estimate.  We  have  found  that  this  procedure 
produces  very  accurate  and  consistent  measurements  of  phase  velocity. 
The  residual  phase  is  then  differentiated  to  obtain  an  improved  group 
velocity  estimate,  and  the  spectral  amplitudes  of  the  Rayleigh  waves 
are  found  by  taking  the  Fourier  transform  of  the  windowed,  compressed 
s i gna I . 

After  applying  this  procedure  to  several  seismograms,  we  check 
the  results  to  make  sure  there  are  no  phase-reversed  seismograms  or 
other  inconsistencies  in  the  dispersion  curves.  We  then  average  the 
phase  and  group  velocities  and  spectral  amplitudes  for  the  selected 
seismograms.  This  averaging  serves  three  functions  ~  it  smooths  the 
data,  reduces  the  chance  of  error  from  inaccurate  times  or  locations, 
and  provides  standard  deviations  on  the  data  which  are  used  as  weights 
by  the  inversion  code. 
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Figure  6.  Block  diagram  of  surface  wave  analysis  procedure. 


Quality  control  during  path  correction  processing  is  very 
important,  and  is  the  most  time  consuming  part  of  the  procedure.  The 
following  procedure  was  used  to  select  the  data  used  for  path 
correction  processing.  First,  the  data  for  each  station  were  plotted 
and  low  signal-to-noise  ratio  data  were  eliminated  (Figure  7).  The 
selected  data  where  replotted  and  marked  at  the  3  km/sec  group 
velocity  point  on  the  seismograms  to  help  identify  timing  errors  in 
the  data  (Figure  8).  This  step  was  essential  for  the  hand-digitized 
WWSSN  data.  Ten  percent  of  the  data  was  found  to  contain  timing 
errors.  Most  of  these  are  minute  errors,  however,  and  for  these 
seismograms  we  simply  adjusted  the  start  time  by  60  seconds.  More 
troublesome  were  smaller  errors  of  10  or  20  seconds.  Since  we  could 
not  easily  guess  the  error  in  these  seismograms,  they  were  not  used  in 
the  data  processing. 

Since  it  was  necessary  to  process  a  large  amount  of  data  (*•  500 
seismograms  for  ~  50  paths),  the  interactive  procedure  described  by 
Stevens,  et  aL.  (1982)  was  partially  automated.  A  trial  set  of  input 
parameters  for  the  data  analysis  code  TELVEL  were  selected.  These 
parameters  are  listed  in  Table  1.  The  variable  parameters  are  the 
frequency  band  and  filter  Q  for  narrow-band  filtering,  the  windowing 
interval  for  phase-matched  filtering,  output  frequencies,  and  2w 
branch  for  phase  velocity.  The  phase  velocity  branch  was  found  by 
estimating  the  phase  velocity  at  one  frequency  and  finding  the  branch 
closest  to  this  value.  The  values  listed  in  Table  1  were  based  on  our 
experience  with  digital  data  (where  the  correct  branch  can  often  be 
found  by  the  long-period  behavior),  on  worldwide  average  values  from 
Oliver  (1962)  and  on  phase  velocities  measured  in  the  western  United 
States  by  McEvilly  (1964).  One  seismogram  for  each  path  was  processed 
interactively  to  make  sure  the  input  parameters  were  appropriate. 

After  the  initial  processing  it  is  necessary  to  select  the  most 
consistent  results  and  average  them  for  inversion.  To  check  this,  we 
plot  all  of  the  phase  and  group  velocity  curves  to  identify 
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Figure  7.  Six  seismograms  at  WWSSN  Station  COL.  Data  processing  starts  by  selecting  data  with  good 
signal/noise  ratio.  The  third  and  fourth  seismograms  (explosions  Cheshire  and  Estuary) 
were  used  in  path  correction  processing. 


Figure  8.  Six  selected  seismograms  for  Station  COL.  Vertical  dashed  line  corresponds  to  a  group- 
velocity  of  3  km/sec.  Notice  the  minute  error  in  seismogram  6  (Estuary). 


TABLE  1 


INITIAL  PARAMETERS  FOR  TELVEL 


Station 

Type 

and 

Distance 

Narrow-Band 

F i i ter i ng 

Phase-Matched 

F i 1  ter i ng 
Window 

Ha  If -Width 
Seconds 

Phase- 
Ve 1 oc i ty 
Branch 

Output 

Frequencies 

FI  F2  q 

T 

F  C 

FI  F2 

SRO 

< 

5000 

km 

0.01 

0.10 

10 

75 

0.015 

4.00 

0.0125 

0.09 

SRO 

> 

5000 

km 

0.01 

0.10 

20 

75 

0.015 

4.00 

0.0125 

0.09 

WWSSN 

< 

1250 

km 

0.01 

0.18 

8 

40 

0.040 

3.65 

0.025 

0.15 

WWSSN 

< 

2750 

km 

0.01 

0.15 

10 

50 

0.030 

3.80 

0.020 

0.125 

WWSSN  >  2750  km  0.01  0.15  12 


60 


0.020  3.95  0.015  0.125 


inconsistent  data.  There  is  almost  always  some  bad  data  resulting 
from  interference,  digitization  errors,  phase  reversals,  incorrect 
phase  velocity  branches,  or  timing  errors.  In  the  upper  figure  in 
Figure  9,  we  show  phase  and  group  velocities  for  all  seismograms 
processed  at  Station  LON.  Three  seismograms  were  rejected  on  the 
basis  of  the  figure.  The  phase  velocity  curve  which  is  offset  from 
the  rest  resulted  from  a  10  second  timing  error.  The  two  curves  which 
are  inconsistent  at  long-periods  are  the  result  of  digitization 
errors.  This  is  a  common  problem  in  hand-digitized  data  when  there  is 
a  strong  Airy  phase  -  narrow-band  filtering  finds  the  Airy  phase  even 
at  the  longest  periods,  so  a  nearly  constant  group  velocity  curve 
results.  We  have  never  seen  this  effect  in  digital  data.  The  lower 
figure  shows  the  average  phase  and  group  velocities  together  with 
_♦  one  standard  deviation  curves  for  the  remaining  nine  seismograms  at 
Station  LON. 

In  Figure  10,  we  show  the  final  phase  and  group  velocities  found 
for  NTS  to  WWSSN  stations,  NTS  to  $R0  stations,  and  East  Kazakh  to  SRO 
stations.  The  phase  and  group  velocities  are  very  consistent  for 
different  events  observed  along  a  given  path.  Typical  phase  velocity 
standard  deviations  are  0.1  to  0.2  percent,  while  group  velocity 
standard  deviations  are  a  few  tenths  of  a  percent.  The  standard 
deviations  of  all  quantities  increase  near  the  ends  of  the  observable 
frequency  band,  and  also  near  spectral  dips. 


18 


Phase  (top)  and  Group  Velocity  Phase  (top)  and  Group  Velocity 

(km/sec)  (km/sec) 


rrj  v«»v 


>  HBHSBHHBinBBaiBnaPBI 


■■■W  ■nmvHDi  npan 


_ i _ i _ i _ i _ i _ i _ 

0.05  0.10 

Frequency  (Hz) 


0.15 


Figure  9.  Phase  and  group  velocities  for  path  from  NTS  to  LON  (1155  km),  top 
figure  shows  velocities  for  each  seismogram.  Bottom  figure  shows 
average  (solid)  and  *  Iff  curves  after  deleting  data  with  timing  and 
digitization  errors.  Origin  id  numbers  are  marked  for  inconsistent 
events . 
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Figure  10. 


Phase  and  group  velocities  for  paths  from  NTS  to  WWSSN  stations 
(top),  NTS  to  SRO  stations  (middle)  and  East  Kazakh  to  SRO  stations 
(bottom)  Solid  lines  are  phase  velocities;  dashed  are  group 
velocities.  The  three  paths  with  high  phase  velocity  in  the  NTS  to 
SRO  figure  are  predominantly  oceanic  paths.  The  three  dispersion 
curves  with  very  low  group  velocities  at  0.03  Hz  in  the  East  Kazakh 
figure  correspond  to  the  paths  with  thick,  low  velocity  crusts  to  the 
south  of  the  test  site. 
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V.  INVERSION  FOR  SHEAR  VELOCITY  STRUCTURE 


The  shear  velocity  structure  is  found  by  simultaneously 
inverting  the  phase  and  group  velocities.  The  method  used  is  similar 
to  the  method  used  by  Bache,  et  al.  (1978a).  The  inversion  program 
attempts  to  fit  the  data  with  a  smooth  model  by  minimizing  the 
integral  of  ld/7/dzl  over  the  structure,  except  across  discontinuities 
at  depths  specified  by  the  user.  The  inversion  method  is  iterative, 
with  the  inversion  model  for  each  iteration  becoming  the  starting 
model  for  the  next  iteration.  The  final  results  are  completely 
independent  of  the  starting  model,  and  the  initial  starting  model  is 
necessary  only  to  allow  the  calculation  of  partial  derivatives  of 
phase  and  group  velocities  with  respect  to  model  parameters. 

Inversion  is  performed  for  shear  velocity  only,  as  there  is  not 
enough  information  in  the  fundamental  mode  phase  and  group  velocities 
to  determine  the  compressions  I  velocities  and  densities  at  the  same 
time.  We  use  data  on  crust  and  upper  mantle  properties  from  Dobrin 
(1976)  to  constrain  these  parameters.  Compressions!  velocity  is 
constrained  to  be  consistent  with  a  Poisson's  ratio  of  0.27,  while 
density  is  constrained  by  a  Birch’s  law  given  by  p  s  0.65/7  ♦  400  (in 
MKS  un i ts) . 

A  20  layered  model  is  used  for  the  inversion.  The  thickness  of 
the  layers  and  depth  of  the  model  are  derived  from  the  frequency 
content  of  the  data.  The  maximum  depth  is  approximately  equal  to  the 
longest  wavelength  c/f,  and  layers  are  thickest  near  the  maximum  depth 
and  become  thinner  near  the  free  surface.  The  inversion  is  damped,  so 
that  the  model  chosen  has  approximately  five  or  six  degrees  of  freedom 
depending  on  the  quality  of  the  data.  This  allows  an  excellent  fit  to 
the  data  in  most  cases,  and  produces  a  smooth  model  free  of  the 
spurious  oscillations  which  would  result  from  allowing  more  degrees  of 
f  reedom . 
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Disconti nuties  in  the  structure  are  explicitly  allowed  at  the 
crust-mantle  boundary  and  in  some  cases  at  shallower  depths  to  allow 
for  near  surface  sedimentary  layers.  The  depths  of  the 
discontinuities  are  found  in  the  following  way.  We  first  perform  an 
inversion  to  find  the  best  model  with  no  discontinuities.  Then  we 
insert  discontinuities  at  each  layer  boundary,  invert,  and  search  for 
minima  in  the  resulting  data  misfit.  Discontinuities  are  inserted  at 
these  depths,  and  a  final  inversion  is  performed  with  the 
d i scont i nu i t i es  i nc I uded . 

In  Figure  11,  we  show  the  data  misfit  minima  found  using  the 
initial  structure  (with  no  discontinuities)  shown  in  Figure  12.  Data 
misfits  are  shown  for  several  degrees  of  freedom.  The  data  misfit  is 
normalized  to  the  standard  deviations  on  the  data,  so  a  data  misfit  of 
1.0  means  that  the  average  misfit  between  the  calculated  and  observed 
phase  and  group  velocities  is  one  standard  deviation.  Minima  are 
found  close  to  7.5  km  and  40  km,  so  discontinuities  are  inserted  at 
these  depths.  The  final  model  is  also  shown  in  Figure  12. 
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Figure  11.  Data  Misfit  minima  for  NTS-AAM.  Discontinuities  were  inserted  after 
each  layer.  Minima  in  data  misfit  correspond  to  discontinuities  in 
structure.  Curves  are  marked  for  different  degrees  of  freedom.  The 
final  structure  (see  next  figure)  has  discontinuities  at  7.5  and  40 


Shear  Velocity  (km/sec) 


VI.  INVERSION  FOR  Q  STRUCTURE  AND  SCALAR  MOMENT 


The  final  step  in  the  inversion  procedure  is  the  inversion  of 
the  average  spectral  amplitude  for  Q  and  moment.  The  average  spectral 
amplitude  is  the  linear  average  of  the  spectra  of  the  best  events 
recorded  for  a  given  path,  defined  by 

i  N 

u(w)  =  jj  ZZMw)  (3) 

"  i=l  ' 


where  u-,  are  the  individual  spectra  and  N  is  the  number  of 
frequencies.  The  relative  moments  of  the  individual  spectra  may  be 
found  using  the  least  squares  relation 
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so  once  the  absolute  moment  corresponding  to  the  average  spectrum  is 
determined,  the  moments  of  all  the  events  are  determined.  Standard 
deviations  on  the  normalized  spectra  (spectra  divided  by  their 
relative  moments)  are  used  as  weights  in  the  inversion.  The  standard 
deviations  are  defined  by 
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In  Figure  13,  we  show  the  spectra  of  individual  seismograms  and  the 
average  spectrum  with  standard  deviations  for  the  path  from  NTS  to 
Station  LON. 
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Figure  13.  Spectral  amplitudes  for  NTS-LON  path.  Top  figure  shows  individual 
spectra  for  nine  explosions.  Bottom  figure  shows  average  spectra  and 
standard  deviation  as  defined  in  the  text. 


Rewriting  Equation  (1),  we  define  the  data  used  for  the 
inversion  to  be  the  ratio  between  the  synthetic  spectral  amplitude 
S (u)  (for  unit  moment,  without  including  attenuation)  and  the  observed 
amplitude  spectrum  u (u) : 


-  In  =  7 (u)  -  -  In  M  '  . 

r  u(«)  ' v  '  r  0 


S(w)  is  defined  by 

S(«)  =  Sj(»)S2(»)/J.e  sin  (r/.e)  .  (7) 

The  attenuation  coefficients  7(0)  can  be  expanded  in  terms  of  the  Q 
structure  of  the  path  using  the  phase  velocity  partial  derivatives 
(Mitchell,  1975;  Anderson,  et  a L.,  1965).  Q  is  assumed  to  be  depth 

dependent,  but  frequency  i ndependent .  Explicitly,  assuming  no  bulk 
losses 
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where  the  index  i  indicates  the  layer  in  the  inversion  model  and  N  is 
the  total  number  of  layers.  a;,  p\ ,  and  Q;  are  the  compress i ona I 
velocity,  shear  velocity  and  Q,  respectively,  in  each  layer. 

In  the  Q/moment  inversion,  the  amplitude  of  the  spectral  ratio 
in  Equation  (6)  determines  the  moment,  while  the  shape  determines  the 
Q  structure.  We  use  a  smoothness  criterion,  as  in  the  inversion  for 
shear  velocity  structure  to  damp  oscillations  in  the  solution.  Since 
there  is  some  correlation  between  low  velocity  zones  and  low  Q  zones 
in  the  earth,  we  invert  for  6/ Q,  and  find  the  smoothest  0/Q  structure 


consistent  with  the  date.  The  Q/moment  inversion  is  less  robust  than 
inversion  for  structure,  since  there  is  an  additional  free  parameter 
(the  moment)  and  since  spectra  are  less  smooth  than  dispersion  curves. 
As  a  result,  the  solution  must  be  more  heavily  damped  by  allowing 
fewer  degrees  of  freedom  in  the  inversion.  The  minimum  number  of 
degrees  of  freedom  is  two  -  one  to  determine  the  moment  and  one  to 
find  the  average  /J/ Q  value  for  the  path.  With  exactly  two  degrees  of 
freedom,  Q  in  each  layer  is  proportional  to  the  shear  velocity,  with 
the  proportional ity  constant  determined  by  the  inversion.  More 
degrees  of  freedom  allow  more  shape  in  the  Q  structure. 

We  impose  one  additional  constraint  on  the  Q  inversion.  Tests 
on  synthetic  seismograms  show  that  while  the  Q  inversion  procedure 
accurately  reproduces  the  Q  structure  of  the  crust,  lower  Q  zones 
beneath  the  crust  are  not  well  resolved,  and  the  Q  at  depth  is  usually 
overestimated.  This  causes  the  moments  measured  at  distant  stations 
to  be  underestimated.  We  can  reduce  this  problem  by  including  an 
estimated  value  for  the  average  mantle  Q  at  depth  as  a  data  point. 
With  this  constraint,  the  Q  inversion  is  stable  and  physically 
reasonable  Q  structures  are  obtained  with  an  inversion  that  allows 
approximately  3.5  degrees  of  freedom  in  the  model. 


VII.  SHEAR  VELOCITY  INVERSION  RESULTS 


Several  representative  shear  velocity  structures  obtained  by 
inversion  are  shown  in  Figure  14.  Three  continental  and  one 
predominantly  oceanic  structures  are  shown  for  NTS.  The  continental 
structures  are  all  quite  similar.  The  path  to  SHA  has  unusually  low 
velocities  in  the  upper  few  kilometers.  Note  the  corresponding  low 
group  velocities  for  this  path  in  Figure  10.  The  path  from  East 
Kazakh  to  KONO  crosses  a  stable,  relatively  fast  shield  structure. 
The  path  to  MAJO  is  similar,  but  somewhat  slower  and  more  complex. 
The  path  to  SHIO  crosses  the  Tibetan  plateau,  and  the  inferred 
structure  has  a  thick,  low  velocity  crust.  The  path  to  ANTO  has  low 
velocities  in  the  upper  part  of  the  crust. 

Our  main  interest  here  is  in  the  effect  of  these  structures  on 
surface  wave  amplitudes  In  Figure  15,  we  show  the  path  amplitude 
functions  Sj  calculated  using  all  of  the  path  structures.  The  elastic 
structure  of  the  source-to-recei ver  path  has  surprisingly  little 
effect  on  the  spectral  amplitudes.  The  maximum  variation  in  log  (S2) 
is  only  0.15  at  20  seconds.  The  paths  with  the  largest  values  of  S2 
are  the  thick  crust,  low  velocity  paths,  while  the  smallest  values 
occur  along  oceanic  paths. 

These  results  mean  that  the  surface  wave  amplitudes  are  not 
changed  very  much  by  the  elastic  properties  of  the  travel  path.  This 
is,  of  course,  only  true  for  spectral  amplitudes,  not  for  time  domain 
amplitudes.  Differences  in  dispersion  can  cause  large  differences  in 
time  domain  surface  wave  amplitudes,  especially  at  distances  less  than 
30  degrees.  On  long  or  complex  paths,  when  a  plane- layered  model  is 
not  a  good  approximation,  interference  will  cause  more  spectral 
amplitude  variation  than  is  predicted  by  the  smooth  functions  in 
Figure  15. 
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Figure  15.  Path  amplitude  function  So  for  all  paths  from  NTS  and  East  Kazakh. 


VIII.  q  CONSTRAINTS  AND  INVERSION  RESULTS 


To  invert  for  the  path  Q  structures  and  the  explosion  moments, 
we  need  the  excitation  function  Sj  for  the  source  region  structures, 
and  an  estimate  of  mantle  Q  at  depth  to  be  used  as  a  constraint. 

The  Soviet  test  site  is  located  in  the  Kazakh  Fold  System  near 
latitude  50  degrees  north  and  longitude  79  degrees  east.  According  to 
Beliayevsky,  et  aL.  (1968),  the  results  of  a  deep  seismic  sounding 
survey  show  that  the  average  crustal  thickness  in  the  area  is  50  km. 
Close  to  the  test  site,  however,  the  crustal  thickness  decreases  to 
about  45  km.  Compressional  velocities  in  the  crust  are  very  similar 
to  velocities  in  the  stable  shield  and  platform  regions  of  the 
Eurasian  continent,  except  for  the  absence  of  low  velocity  surface 
sediments  in  the  Kazakh  Fold  System.  Grand,  et  aL.  (1984)  show  that 
the  shear  velocity  in  the  crust  near  the  source  region  is  relatively 
slow,  similar  to  the  shear  velocity  of  the  platforms  to  the  south  and 
east,  but  slower  than  the  velocity  along  the  shield  paths  to  Europe. 

For  our  model  of  the  East  Kazakh  source  region,  we  use  a 
modification  of  the  structure  inferred  for  the  path  to  MAIO.  This 
path  crosses  the  Turanian  Plateau  southwest  of  the  test  site.  The 
crustal  thickness  of  this  structure  is  45  km.  We  have  increased  the 
inferred  velocities  in  the  top  15  km  (see  Figure  16)  to  remove  the 
influence  of  the  low  velocity  sediments  in  the  plateau.  The  structure 
in  the  upper  crust  is  then  consistent  with  velocities  from 
Beliayevsky,  et  aL.  (1968).  We  have  reduced  the  velocity  in  the  top 
1.5  km  to  5.0  km/sec  (compressional)  and  2.7  km/sec  (shear)  with  a 

O 

density  of  2.7  gm/cm  .  These  values  are  appropriate  for  weathered 
hard  rock  such  as  granite. 
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The  models  used  for  NTS  and  East  Kazakh  are  listed  in  Tables  2 
and  3  and  shown  in  Figure  17.  The  NTS  model  is  based  on  the  NTS-TUC 
crustal  structure  of  Bache,  et  aL.  (1978a),  merged  with  the  upper 
mantle  structure  of  Anderson  and  Hart  (1976,  1978).  The  shallow 
structure  (upper  4  km)  is  the  Pahute  Mesa  structure  from  Bache,  et  aL. 
(1975).  The  shear  velocity  at  the  source  depth  of  1  km  is  2000  m/sec, 
and  the  shear  velocity  in  the  top  500  m  is  1000  m/sec. 

The  excitation  functions  for  the  two  source  regions  are  shown  in 
Figure  18.  The  excitation  function  for  NTS  is  greater  than  the 
excitation  function  for  East  Kazakh  at  all  frequencies  of  interest, 
and  particularly  at  higher  frequencies.  We  therefore  expect  an 
explosion  at  NTS  to  generate  larger  surface  waves  than  an  explosion  at 
East  Kazakh  with  the  same  moment. 

We  also  need  an  estimate  of  worldwide  average  mantle  Q  to  use  as 
a  constraint  on  the  inversions,  so  we  have  searched  through  the 
literature  to  find  a  representati ve  set  of  Q  estimates.  Anderson  and 
Hart  (1978)  find  a  worldwide  average  Q  of  90  to  110  between  80  and  300 
km  depth.  Canas  and  Mitchell  (1978)  find  a  minimum  mantle  Q  of  100  to 
200  in  the  Pacific  Ocean  with  the  minimum  at  about  150  km  depth,  while 
Canas  and  Mitchell  (1981)  find  a  minimum  mantle  Q  of  60  to  150  in  the 
Atlantic  Ocean  at  100  to  150  km  depth.  Cara  (1981)  finds  a  mantle  Q 
of  80  to  100  between  100  and  300  km  depth  for  paths  across  the  Pacific 
Ocean.  Lee  and  Solomon  (1979)  find  very  low  mantle  Q  values  of  about 
20  in  the  Western  United  States  between  80  and  160  km,  25  in  the  east- 
central  United  States  below  135  km,  and  70  in  the  Central  Pacific 
below  60  km.  Singh  (1982)  finds  a  mantle  Q  of  about  50  in  southern 
Asia.  Patton  (1980)  finds  mantle  Q  values  of  60  to  150  in  western 
Eurasia  at  depths  below  50  km,  100  to  150  in  eastern  Eurasia  at  depths 
below  100  km,  and  300  to  500  in  the  stable  platforms  of  Eurasia  at 
depths  below  50  km.  There  is  a  considerable  amount  of  uncertainty  in 
all  of  these  results,  and  considerable  variation  in  the  estimates. 
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TABLE  2 


Thickness 

(k«) 


NTS  SOURCE  REGION  STRUCTURE 


P-V®loc»ty  S-Velocity 

(km/sec)  (km/sec) 


Shear  Velocity  (km/cec)  Shear  Velocity  (knt/eec) 


In  the  inversions  described  in  this  paper,  we  set  0/Q  equal  to 
40  with  a  standard  deviation  of  1  between  120  and  150  km  depth.  Since 
the  shear  velocity  at  this  depth  is  between  4000  and  5000  m/sec,  Q  is 
constrained  to  be  slightly  more  than  100  in  this  depth  range.  The 
inversion  algorithm  treats  this  "constraint*  as  a  data  point  and  tries 
to  find  a  Q  structure  that  both  fits  the  observed  spectral  amplitude 
data  and  has  the  preferred  Q  value  in  the  specified  depth  range.  The 
Q  value  from  the  inversion  may  vary  from  the  preferred  value  if  the 


data  requires  it. 

Several  representative  Q  structures  are  shown  in  Figure  19. 
Inferred  crustal  Q  is  quite  variable,  with  very  low  values  for  some 
complex  paths  such  as  East  Kazakh  to  ANTO  in  Turkey.  Seismograms 
recorded  at  KQNO  (from  East  Kazakh)  are  rich  in  high  frequencies, 
which  may  imply  a  high  crustal  Q  along  this  path. 

In  Figure  20,  we  show  the  inferred  attenuation  coefficients  for 
all  paths.  In  Table  4,  we  have  listed  the  average  20  second  values  of 
attenuation  coefficients,  together  with  the  structural  parameters,  for 
paths  from  NTS  to  WWSSN  stations  and  SRO  stations,  for  the  three 
oceanic  paths  out  of  NTS,  for  paths  from  East  Kazakh  to  SRO  stations 
and  for  all  paths.  Differences  in  attenuation  may  be  responsible  for 
large  amplitude  variations.  Rayleigh  waves  observed  along  paths  with 
attenuation  coefficients  of  1.0  x  1.0'^/meters  and  3.0  x  10~^/meters 
(the  *  one  standard  deviation  values  for  all  paths)  will  differ  by 
0.09  magnitude  units  at  1000  km,  0.43  magnitude  units  at  5000  km,  and 
0.87  magnitude  units  at  10,000  km. 


Q  Structures 


structures  found  by  Q/moment  inversion  for  selected  paths. 


Attenuation  Coefficients  (y)-AII  Paths 


Figure  20.  Attenuation  coefficients  for  all  paths  from  NTS  and  East  Kazakh. 


TABLE  4 

AVERAGE  20  SECOND  VALUES  FOR  NTS  AND 
EAST  KAZAKH  PATHS  (MKS  UNITS) 


IX.  ESTIMATED  MOMENTS  FOR  NTS  AND  EAST  KAZAKH  EXPLOSIONS 

We  have  estimated  explosion  moments  for  the  40  NTS  explosions 
and  16  East  Kazakh  explosions  which  were  recorded  at  two  or  more 
stations.  The  estimated  moments  (in  MKS  units,  newton-meters), 
together  with  their  network  standard  deviations  and  the  number  of 
stations  used  in  the  average  are  given  in  Tables  5  and  6.  The  yields 
given  m  Table  5  are  from  Marshall,  et  aL.  (1979).  The  network 
standard  deviations  in  (log)  moment  are  quite  small,  about  0.1.  Even 
for  recent  NTS  explosions,  which  include  data  from  several  distant  SRO 
stations,  the  standard  deviations  are  only  about  0.15. 

As  a  check  on  the  results,  we  calculated  residual  station 
corrections  for  all  paths.  These  are  listed  in  Tables  7  and  8.  These 
station  corrections  are  designed  to  minimize  the  variance  in  moment 
over  the  network  at  stations  observing  each  event.  Most  of  the 
station  corrections  are  small,  particularly  for  shorter  paths. 
Eighteen  of  the  24  WWSSN  stations  have  station  residuals  of  0.1  or 
less.  The  station  corrections  increase  with  distance  from  the  source, 
which  is  to  be  expected  because  of  the  increasing  complexity  of  the 
travel  paths,  but  the  largest  corrections  for  some  very  complex  paths 
are  only  about  0.3. 

Sykes  and  Cifuentes  (1984)  measured  Ms  values  for  a  similar  set 
of  East  Kazakh  explosions.  Our  moment  estimates  are  very  consistent 
with  their  reported  Ms  values.  Using  the  ten  normal  events  that  the 
two  data  sets  have  in  common  we  find  the  following  moment/Ms  relation: 

log  M0  =  Ms  ♦  11.86  *  0.03 

where  the  *  0.03  is  the  standard  deviation  over  all  ten  events. 

It  is  interesting  to  compare  the  results  for  recent  events  with 
theoretical  explosion  models.  Since  April  1976,  underground 
explosions  have  been  limited  to  a  yield  of  150  kilotons  by  the 
Threshold  Test  Ban  Treaty.  In  Figure  21,  we  show  the  estimated 


TABLE  5 


NTS  EXPLOSIONS 


Number 


Event 

Date 

of 

Area*  Stations 

Log 

Moment 

Standard 

Deviation 

Corduroy 

12/03/65 

Y 

5 

16.10 

0.14 

Dumont 

05/19/66 

Y 

5 

15.98 

0.13 

Pi ledr i ver 

06/02/66 

C 

6 

15.84 

0.12 

Commodore 

05/20/67 

Y 

5 

16.27 

0.26 

Scotch 

05/23/67 

P 

7 

16.46 

0.18 

Zaza 

09/27/67 

Y 

6 

16.32 

0.14 

Lanpher 

10/18/67 

Y 

4 

15.87 

0.06 

St i nger 

03/22/68 

P 

3 

16.17 

0.09 

Sled 

08/29/68 

P 

10 

16.22 

0.09 

Nogg i n 

09/06/68 

Y 

6 

16.05 

0.06 

Purse 

05/07/69 

P 

9 

16.42 

0.06 

Ti jeras 

10/14/70 

Y 

5 

16.01 

0.13 

Oscuro 

09/21/72 

Y 

3 

15.99 

0.11 

Almendro 

06/06/73 

P 

10 

16.78 

0.11 

Escabosa 

07/10/74 

Y 

7 

16.14 

0.15 

Portmanteau 

08/30/74 

Y 

2 

16.08 

0.00 

Topgal lant 

02/28/75 

Y 

2 

15.83 

0.09 

Tybo 

05/14/75 

P 

10 

16.50 

0.10 

Sti (ton 

06/03/75 

P 

5 

16.17 

0.12 

Mi zzen 

06/03/75 

Y 

5 

16.02 

0.15 

Mast 

06/19/75 

P 

12 

16.58 

0.15 

Cheshire 

02/14/76 

P 

4 

16.79 

0.39 

Estuary 

03/09/76 

P 

4 

16.79 

0.05 

Poo  1 

03/17/76 

P 

3 

16.52 

0.09 

Strait 

03/17/76 

Y 

3 

16.28 

0.15 

Marsi 1 ly 

04/05/77 

Y 

2 

16.01 

0.06 

Sandreef 

11/09/77 

Y 

15 

16.27 

0.15 

Fara 1 1  ones 

12/14/77 

Y 

2 

16.22 

0.02 

Iceburg 

03/23/78 

Y 

2 

16.15 

0.07 

Lowba 1 1 

07/12/78 

Y 

4 

15.96 

0.07 

Pan  i  r 

08/31/78 

P 

8 

15.97 

0.15 

Rummy 

09/27/78 

Y 

16 

16.12 

0.15 

Farm 

12/16/78 

P 

7 

15.89 

0.18 

Pepato 

06/11/79 

P 

17 

16.07 

0.13 

Hearts 

09/06/79 

Y 

16 

16.12 

0.15 

Sheepshead 

09/26/79 

P 

14 

15.95 

0.14 

Pyramid 

04/16/80 

Y 

5 

15.65 

0.08 

Colwick 

04/26/80 

P 

7 

15.91 

0.18 

Kash 

06/12/80 

P 

8 

16.16 

0.09 

Taf  i 

07/25/80 

P 

7 

16.10 

0.09 

*Y  =  Yucca 

Flat,  P  =  Pahute  Mesa, 

C  =  Cl 

max  Stock 

Remarks 


Yield  62  KT,  Reversed 
Yield  250  KT 
Yield  155  KT 


150  KT  Limit  Begins 
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TABLE  6 

EAST  KAZAKH  EXPLOSIONS 

Number 

of 

Log  Standard 

Stations 

Moment  Dev i at i on 

Remarks 

07/04/76 

3 

15.73 

0.18 

06/11/78 

5 

15.94 

0.04 

:  09/15/78 

6 

15.74 

0.10 

11/29/78 

9 

15.75 

0.11 

06/23/79 

10 

15.83 

0.09 

08/04/79 

10 

15.94 

0.12 

08/18/79 

9 

15.69 

0.20 

Reversed 

10/28/79 

6 

15.85 

0.08 

i  12/02/79 

9 

15.96 

0.11 

'  12/23/79 

3 

15.66 

0.20 

Large  body  waves 

09/14/80 

10 

15.86 

0.09 

Reversed 

10/12/80 

7 

15.94 

0.09 

.  04/22/81 

7 

15.92 

0.09 

•  09/13/81 

4 

15.87 

0.32 

Very  low  at  KAAO 

10/18/81 

3 

15.91 

0.13 

12/27/81 

5 

15.92 

0.09 

04/25/82 

5 

15.85 

0.10 

j  12/05/82 

5 

15.92 

0.12 

moments 
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MARSILLY 

5.6 

SANOREEF 

5.7 

FARALLONES 

5.7 

ICEBERG 

5.6 

LOWBALL 

5.5 

PANIR 

5.6 

RUMMY 

5.7 

FARM 

5.5 

PEPATO 

5.5 

HEARTS 

5.8 

SHEEPSHEAD  5.6 

PYRAMID 

5.3 

COLWICK 

5.4 

KASH 

56 

TAFI 

5.5 

9/15/78 

6.< 

11/29/78 

6.i 

8/4/79 

6. 

8/18/79 

6. 

10/28/79 


12/2/79 

6. 

12/23/79 

6 

9/14/80 

6.2 

10/12/80 

5.9 

4/22/81 

5.9 

9/13/81 

6.0 

12/27/81 

6.2 
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moments  for  the  largest  explosions  at  NTS  and  East  Kazakh  since  this 
date.  Using  a  Muel ler-Murphy  explosion  model  (Mueller  and  Murphy, 
1971),  the  log  moment  of  a  150  kiloton  explosion  in  a  tuff/rhyolite 
medium  should  be  16.24.  This  agrees  very  well  with  the  upper  limit  of 
the  recent  NTS  explosion  estimated  moments. 


A  surprising  result  is  the  low  values  of  moment  estimated  for 
the  East  Kazakh  explosions.  None  of  the  explosions  has  an  estimated 
log  moment  greater  than  16.0.  The  upper  limit  of  the  log  moments  for 
East  Kazakh  is  about  0.2  less  than  the  upper  limit  of  log  moments  for 
NTS.  This  contiasts  dramatically  with  the  NEIS  m^  values  for  the 
explosions  (also  shown  in  Figure  21).  Several  of  the  East  Kazakh 
explosions  have  m^  greater  than  6.1,  while  none  of  the  NTS  explosions 
has  an  m^  greater  than  5.8.  Comparison  with  theoretical  explosion 
models  makes  the  contrast  even  larger.  It  is  commonly  believed  that 
the  source  medium  at  the  East  Kazakh  test  site  is  a  hard  "granite- 
like”  material  (e.g.,  Alewine  and  Bache,  1983),  but  the  log  moment  for 
a  150  kiloton  Muel ler-Murphy  explosion  in  granite  is  16.52,  more  than 
0.5  greater  than  the  largest  estimated  moment  at  East  Kazakh. 
Furthermore,  the  difference  between  m^  and  M$  due  to  local  elastic 
properties  at  the  source  should  be  larger  for  a  low  velocity  material 
then  for  a  high  velocity  material  (Stevens  and  Day,  1985). 

Several  factors,  including  tectonic  strain  release,  differences 
in  body  wave  attenuation,  or  an  unusual  source  medium  could  account 
for  these  differences.  Correction  for  tectonic  strain  release  could 
increase  the  moments  of  some  events  to  the  level  of  the  largest  NTS 
explosions  (Stevens,  1982;  Sykes  and  Cifuentes,  1984),  but  could  not 
increase  them  to  the  level  expected  for  a  150  kiloton  granite  source. 
The  East  Kazakh  explosions  are  therefore  not  consistent  with  a 
"granite- 1  ike”  source  medium.  One  possible  explanation  for  the  large 
body  waves  and  small  surface  waves  from  East  Kazakh  explosions  would 
be  the  presence  of  a  source  medium  at  the  test  site  that  resulted  in  a 
highly  peaked  source  spectrum  with  a  large  amplitude  near  one  Hz  and 
an  unusually  small  long  period  level. 
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X.  THE  EFFECT  OF  SOURCE  REGION  STRUCTURE  ON  MOMENT  ESTIMATES 


It  is  important  to  assess  the  consequences  of  the  constraints 
and  assumptions  that  were  used  in  making  our  moment  estimates.  The 
two  most  important  constraints  were  the  model  for  the  elastic 
structure  of  each  source  region  and  the  constraint  placed  on  the  upper 
mantle  Q.  In  this  section  and  the  following  section,  we  examine  the 
variations  in  moment  estimates  that  result  from  changes  in  these 
assumptions. 


The  East  Kazakh  source  region  structure  is  based  on  the  best 
available  information  about  the  region,  but  this  information  cannot  be 
regarded  as  definitive  so  we  want  to  estimate  the  maximum  reasonable 
error  in  our  excitation  function  for  this  region.  In  Figure  22,  we 
show  the  excitation  functions  calculated  from  all  of  the  earth 
structures  that  we  have  found  by  surface  wave  inversion.  There  is  a 
substantial  amount  of  variation  in  the  excitation  functions  with  a 
total  variation  of  about  a  factor  of  two.  The  smallest  excitation 
functions,  however,  belong  to  oceanic  and  mixed  paths,  and  the  largest 
are  for  thick  crust,  low  velocity  paths  to  the  southeast  of  East 
Kazakh.  In  Figure  23  we  show  the  excitation  functions  for  all  of  the 
continental  paths  out  of  East  Kazakh.  The  variation  in  is  reduced 
substantially.  The  lowest  amplitude  excitation  function  belongs  to 
the  high  velocity  shield  path  from  East  Kazakh  to  station  KONO. 


To  estimate  the  uncertainty  in  the  excitation  function,  we  used 
the  East  Kazakh  to  KONO  structure,  with  the  surface  velocities 
modified  to  be  the  same  as  in  the  structure  used  for  our  Q/moment 
inversions  described  earlier.  We  calculated  the  excitation  function 
for  this  structure  and  inverted  the  spectra  for  all  paths  out  of  East 
Kazakh  to  obtain  new  moment  estimates.  The  log  moments  were  increased 
by  0.04  on  the  average. 


Source  Region  Excitation  Function  (SI )- All  Paths 


Figure  22.  Excitation  function  Sj  calculated  using  the  structures  found  for  all  paths  from  NTS  and 
East  Kazakh. 


Figure  23.  Excitation  functions  calculated 


The  most  important  feature  of  the  NTS  source  region  structure 
is  the  shallow  shear  velocity.  In  Figure  24  we  show  the  variation  of 
the  source  region  excitation  function  calculated  at  a  frequency  of 
0.05  Hz  as  a  function  of  the  shear  velocity  at  the  source  depth  of  one 
kilometer.  In  these  calculations,  the  shear  velocity  in  the  upper  1.5 
kilometers  of  the  NTS  source  region  model  was  varied,  the 
compress i ona I  velocity  was  set  using  a  constant  Poisson’s  ratio  of 
0.27  and  the  density  in  this  layer  was  held  fixed.  The  phase  and 
group  velocities  at  this  low  frequency  depend  only  weakly  on  the  shear 
velocity,  but  the. excitation  function  S^  is  quite  sensitive  to  the 
shallow  velocity  structure,  varying  approximately  as  S^  ~ 

The  shear  velocity  at  the  source  depth  used  in  the  earlier 
moment  estimates  for  NTS  was  2.0  km/sec.  Bache,  et  aL.  (1975)  give 
shear  velocities  in  tuff  and  rhyolite  at  NTS  ranging  from  1.5  to  2.5 
km/sec,  so  this  is  a  good  average  value  for  the  shear  velocity  at  NTS. 
The  lowest  hard  rock  velocities  at  NTS  are  found  in  dry  tuff,  which 
have  shear  velocities  of  about  1.0  km/sec.  To  put  an  upper  bound  on 
the  NTS  source  region  excitation  function,  we  modified  the  NTS  model 
to  have  a  shear  velocity  in  the  upper  1.5  km  of  1.0  km/sec,  calculated 
the  excitation  function  for  this  structure,  and  reinverted  all  of  the 
spectra  for  NTS  paths  to  obtain  new  moment  estimates.  Log  moment 
estimates  were  reduced  by  0.08  on  the  average. 

The  excitation  functions  for  NTS  and  East  Kazakh  are  shown  in 
Figure  25  together  with  the  excitation  functions  calculated  from  the 
two  perturbed  structures  described  here.  The  excitation  functions  are 
very  frequency  dependent,  so  it  is  not  possible  to  tell  at  a  glance 
what  the  effect  of  the  excitation  function  is  on  the  moment  estimates. 
To  estimate  the  effect  of  the  different  source  regions  used  for  NTS 
and  East  Kazakh,  we  re  inverted  the  spectra  for  East  Kazakh  paths  using 
the  NTS  source  region  model .  The  log  moments  were  decreased  on  the 
average  by  0.02. 


20  Second  Sensitivity  to  Shallow  Shear  Velocity 


Figure  24.  The  dependence  of  the  20  second  excitation  function 


Low  Velocity  NTS 


Excitation  functions  calculated  from  NTS  and  East  Kazakh  source  region  structures.  The 
lower  dashed  line  is  the  excitation  function  derived  from  the  East  Kazakh  to  KONO  path 
structure.  The  upper  dashed  line  is  the  excitation  function  for  the  NTS  structure  with 
the  shear  velocity  at  the  source  depth  reduced  to  1.0  km/sec.  Estimates  for  Alog  mb  were 
obtained  by  inverting  spectra  from  NTS  and  East  Kazakh  paths  for  Q  and  moment  using  each 


Correction  for  source  region  structure  reduces  estimated 
moments  for  NTS  explosions  by  approximately  0.02  relative  to  estimated 
moments  for  East  Kazakh  explosions.  This  correction  would  be  larger 
if  either  the  NTS  near  surface  velocities  were  lower  than  the  values 
used,  or  if  the  East  Kazakh  crustal  velocities  were  higher  than  those 
used  in  the  model.  Note  also  that  the  correction  to  the  spectral 
amplitudes  increases  with  frequency.  The  correction  to  a  20  second 
spectral  magnitude  should  be  approximately  0.09.  The  correction  to 
moment  is  smaller  because  the  moment  estimate  depends  more  strongly  on 
long  period  amplitudes  than  on  short  period  amplitudes,  and  because 
there  is  a  tradeoff  between  the  source  region  correction  and  the 
attenuation  correction.  For  a  given  observed  Rayleigh  wave  spectrum, 
moment  estimation  using  the  NTS  structure  will  result  in  a  larger 
attenuation  correction  than  moment  estimation  using  the  East  Kazakh 
structure,  since  the  inversion  will  assume  (properly)  that  if  the 
observed  spectrum  had  contained  the  amount  of  high  frequency  energy 
shown  in  the  NTS  excitation  function  in  Figure  25  at  the  source ,  it 
must  have  been  attenuated  by  the  time  it  reached  the  receiver. 

One  final  point  should  be  made  about  the  source  region 
correction.  In  Section  II,  we  defined  S^  so  that  we  find  the 
normalized  moment  Mq  rather  than  the  moment  Mq.  This  choice  is 
appropriate  since  Mq  is  a  representation  of  the  source  of  the  surface 
waves  that  is  independent  of  local  material  properties,  and  since 
there  is  no  evidence  that  Mq  correlates  with  yield  better  than  Mq.  We 
can,  of  course,  convert  Mq  to  Mq  if  local  velocities  are  known  through 
the  relation 


3  P 


2*0 


Since  Poisson's  ratio  is  probably  lower  at  East  Kazakh  than  at  NTS 
(because  of  stronger,  more  competent  rock),  this  conversion  would 
increase  estimated  moments  at  NTS  relative  to  estimated  moments  at 
East  Kazakh. 
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XI.  THE  EFFECT  OF  MANTLE  Q  CONSTRAINTS  ON  MOMENT  ESTIMATES 


A  good  check  on  the  validity  of  the  mantle  Q  constraint  is  the 
variation  of  the  moment  estimates  as  a  function  of  distance.  In 
Figure  26,  we  show  the  station  residuals  for  all  NTS  and  East  Kazakh 
paths  plotted  versus  distance.  The  station  residuals  are  small  at 
close  range,  approximately  0.1  for  most  stations  less  than  5000 
kilometers  from  a  test  site.  The  magnitude  of  the  residuals  increases 
with  distance,  as  is  to  be  expected  with  the  increasing  complexity  of 
the  path  at  large  distances,  but  the  residuals  scatter  evenly  around 
zero.  The  lack  of  any  trend  in  the  data  suggests  that  the  attenuation 
correction  has  been  properly  estimated  with  the  constraint  of  /J/ Q  =  40 
between  120  and  150  kilometers  depth. 

In  Figures  27,  28,  and  29  we  show  station  residuals  with 
constraints  of  /J/Q  =25,  50  and  100  respectively.  With  the  lower 
attenuation  of  /7/Q  =25  (Q  “  180),  the  station  residuals  clearly 
increase  with  distance  indicating  that  the  moments  are  undercorrected 
for  attenuation.  The  relative  moments  at  the  two  test  sites  are  not 
changed  significantly,  however.  The  log  moments  are  decreased  by  0.05 
at  NTS  and  by  0.07  at  East  Kazakh.  Increasing  attenuation  slightly 
with  /7/Q  =  50  (Q  ~  90),  the  station  residuals  show  little  trend  with 
distance  and  the  moments  are  increased  by  0.03  at  East  Kazakh  and  by 
0.04  at  NTS.  If  attenuation  is  increased  to  /J/Q  =  100  (Q  ~  45),  then 
the  station  residuals  show  a  strongly  decreasing  trend  with  distance 
indicating  that  the  moments  have  been  overcorrected  for  attenuation. 
In  this  case  the  moments  are  increased  by  0.22  at  NTS  and  by  0.23  at 
East  Kazakh. 

Based  on  these  tests,  we  conclude  that  a  constraint  of 
/?/Q~  40-50  in  the  upper  mantle  is  appropriate.  The  relative  moments 
between  the  two  test  sites  is  insensitive  to  the  choice  of  constraint 
as  long  as  the  same  constraint  is  used  for  both  test  sites  and  the 
distribution  of  stations  is  comparable  for  both  test  sites. 


NTS 


Figure  26.  Station  residuals  with  a  mantle  Q  constraint  of  /7/Q  =  40. 


Station  Residual  Bata/Q  =  25 


with  mantle  Q  constraint  of  /J/Q  =  25. 


Figure  28.  Station  residuals  with  a  mantle  Q  constraint  of  /?/ Q  =  50. 


Station  Residual  Beta/Q  =100 


Figure  29.  Station  residuals  with  a  mantle  Q  constraint  of  p/ Q  =  100. 


throughout  the  earth.  Because  of  the  finite  frequency  range  of  the 
data,  we  cannot  directly  determine  the  mantle  Q,  and  our  approach  has 
been  to  use  an  average  value  as  a  constraint  and  to  accept  some 
scatter  in  the  results,  rather  than  to  try  to  make  uncertain 
corrections  for  variations  in  mantle  Q.  This  is  particularly 
appropriate  for  long  paths,  where  the  Q  co’rrection  is  most  important, 
since  many  different  regions  of  the  earth  are  sampled  by  the  Rayleigh 
wave  between  the  source  and  receiver. 

An  interesting  question  is  the  effect  of  a  small  region  of  the 
travel  path  with  a  mantle  Q  significantly  different  from  the  average 
value.  As  we  mentioned  earlier,  for  example,  interstation  Q  studies 
have  found  much  lower  mantle  Q  values  near  NTS  (e.g.,  Lee  and  Solomon, 
1979)  than  near  East  Kazakh  (e.g.,  Patton,  1980).  If  the  mantle  Q  is 
lower  than  our  average  value,  then  the  attenuation  coefficients  will 
be  underestimated.  From  Equation  (6),  we  find  for  a  moment  estimate 
Mq  made  from  an  observed  spectrum  u(v), 

AlnM0  *  A(r7(«))  (9) 

so  if  7  is  underestimated,  Mq  will  also  be  underestimated  by  an  amount 
that  depends  on  the  length  of  the  travel  path  with  the  different  Q 
value.  Comparing  results  from  our  test  cases  above,  we  find  that 
changing  mantle  /J/ Q  from  40  to  25  decreases  7  by  approximately 
0.3  x  10“4/km  over  the  frequency  range  0.02  to  0.05  Hz.  Changing  /J/Q 
from  40  to  100  increases  7  by  approximately  1.0  x  10“^/km  over  the 
same  frequency  range.  If  l/Q  =  100  were  a  better  value  for  mantle  Q 
for  a  region  approximately  1000  km  around  NTS,  then  our  log  Mq 
estimates  for  NTS  explosions  would  be  low  by  about  0.04.  If  /?/Q  =  25 
were  a  better  value  for  mantle  Q  for  a  region  1000  km  around  East 
Kazakh,  then  our  log  Mq  estimates  would  be  high  by  approximately  0.01. 
Making  this  additional  correction  would  therefore  increase  the  NTS 
moments  relative  to  the  East  Kazakh  moments  by  approximately  0.05. 
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SUMMARY  AND  CONCLUSIONS 


In  summary,  the  most  important  conclusions  of  this  report  are: 

*  The  path  correction  procedure  does  provide  stable  and 
accurate  measurements  of  surface  waves. 

*  Typical  network  scatter  in  log  moment  is  about  0.1 
compared  to  0.2  or  more  for  network  scatter  in  Ms. 

*  Dispersion  is  not  a  problem  since  this  is  a  spectral 

method.  Data  from  stations  at  close  range  (less  than 

30  degrees)  can  be  used  without  difficulty. 

*  The  spectrum  is  flattened  by  application  of  the  path 

correction,  so  the  results  are  nearly  frequency 
independent.  Moments  are  estimated  from  an  average 
over  a  range  of  frequencies. 

*  Little  variability  in  spectral  amplitudes  results  from  the 
elastic  structure  of  the  travel  path. 

*  The  standard  deviation  on  the  logarithm  of  the  path 

amplitude  function  over  all  of  the  paths  processed  is 
only  0.025.  The  extreme  values  differ  by  only  0.15. 

*  This  result  contrasts  dramatically  with  time  domain 
amplitudes  which  have  large  variations  due  to  path 
dependent  dispersion. 

*  The  elastic  structure  of  the  source  region  can 
significantly  affect  surface  wave  amplitudes. 

*  Assuming  that  the  surface  velocities  are  higher  at  East 
Kazakh  than  at  NTS,  surface  wave  amplitudes  should  be 
larger  for  NTS  explosions  than  for  East  Kazakh 
explosions  with  the  same  moment. 
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*  Correction  for  source  region  structure  based  on  the 
best  models  available  for  the  two  test  sites  reduces 
the  estimated  log  moment  of  NTS  explosions  by  0.02 
relative  to  log  moments  for  explosions  at  the  East 
Kazakh  test  site. 

Attenuation  differences  can  cause  large  variations  in 
surface  wave  amplitudes  at  distant  stations.  An  important 
check  on  the  attenuation  correction  is  the  variation  of 
moment  estimates  with  distance. 

*  When  the  Q  inversion  is  performed  with  the  mantle  Q 
constrained  to  a  value  of  approximately  100  at  depths 
of  120  to  150  kilometers,  the  moment  estimates  remain 
constant  (on  the  average)  with  distance. 

*  When  the  Q  inversion  is  performed  with  a  Q  constraint  a 
factor  of  two  higher  or  lower,  the  moment  estimates 
show  a  significant  trend  with  distance  indicating  under 
or  overcorrection  for  attenuation. 

*  The  relative  moments  between  the  two  test  sites  are 
insensitive  to  changes  in  the  Q  constraint. 
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